import cplex
import random
import numpy as np
import time
from pylab import *
from stoch_trnsport_gen import *
from cplex.exceptions import CplexError
import mpl_toolkits.mplot3d.axes3d as p3

NI = 50
NJ = 50
NW = 500
I = range(NI)
J = range(NJ)
O = range(NW)

prob,demand = gen_stoch_trnsport(NI, NJ, NW, eps=0.2)
#prob.solve()

bigm = 0
bigm_w = [0.2 for w in O]


#true_opt_val = prob.solution.get_objective_value()

#print prob.solution.get_objective_value(),
#print prob.solution.progress.get_num_iterations()

bigm_vec = [0.005 * i for i in range(1,21)]

X_bigm_vec = [0.02 * i for i in range(1,21)]
Y_bigm_vec = [0.02 * i for i in range(1,21)]
X_bigm_arr = np.array(X_bigm_vec)
Y_bigm_arr = np.array(Y_bigm_vec)

vals_arr = np.zeros((len(X_bigm_vec), len(X_bigm_vec)))
time_arr = np.zeros((len(X_bigm_vec), len(X_bigm_vec)))

x_bigm_star = 0
y_bigm_star = 0

Z = vals_arr
T = time_arr

for k in range(len(bigm_vec)):

    row_name = []
    col_name = []
    bigm_coefs = []
    rhs_vals = []

    for w in O:
        for j in J:
            row_name.append('scenario' + str(w) + '_' + 'customer' + str(j))
            col_name.append('y' + '_' + str(w))
            bigm_coefs.append( -bigm_vec[k] * demand[w][j] )
            rhs_vals.append( demand[w][j] * (1 - bigm_vec[k]) )

    prob.linear_constraints.set_coefficients( zip(row_name, col_name, bigm_coefs) )
    prob.linear_constraints.set_rhs( zip(row_name, rhs_vals) )
#    prob.linear_constraints.set_coefficients(row_name, col_name, -bigm_vec[k] * demand[w][j])
#    prob.linear_constraints.set_rhs(row_name, demand[w][j] * (1 - bigm_vec[k]) )


    start_time = time.time()
    prob.solve()
    use_time = time.time() - start_time


    print prob.solution.get_objective_value(), 
    print prob.solution.progress.get_num_iterations(),
    print use_time,
    print bigm_vec[k]

    
#vals_arr[l] = prob.solution.get_objective_value()
#time_arr[l, k] = use_time
#vals_arr[l, k] = prob.solution.get_objective_value()
"""
temp_opt_val = prob.solution.get_objective_value()

if abs(temp_opt_val - true_opt_val) <= 0.000001 and use_time <= best_use_time:

    best_use_time = use_time
    x_bigm_star = X_bigm_vec[l]
    y_bigm_star = Y_bigm_vec[k]
""" 

#print X_bigm_vec[l], Y_bigm_vec[k]


